/**
 * 
 */
package ua.lviv.franko.solvers.old.twoDimention;

import java.io.FileInputStream;
import java.io.FileNotFoundException;
import java.util.ArrayList;
import java.util.Iterator;
import java.util.Scanner;

import org.jzy3d.maths.Range;

import ua.lviv.franko.integration.ConstantFunction2D;
import ua.lviv.franko.integration.FuncWithJacobian;
import ua.lviv.franko.integration.Gauss2D;
import ua.lviv.franko.integration.IFunction2D;
import ua.lviv.franko.matrix.RowType;
import ua.lviv.franko.matrix.hash.lu.HashMatrixLUStyle;
import ua.lviv.franko.matrix.hash.lu.LUFactorization;
import ua.lviv.franko.mlsFunctions.oneDimension.weightFunctions.CubicSpline;
import ua.lviv.franko.mlsFunctions.twoDimension.BaseFunction;
import ua.lviv.franko.mlsFunctions.twoDimension.Particle;
import ua.lviv.franko.mlsFunctions.twoDimension.weightFunctions.IWeightFunction2D;
import ua.lviv.franko.mlsFunctions.twoDimension.weightFunctions.SquareSupport;
import ua.lviv.franko.solvers.twoDimention.EquationParameters;
import ua.lviv.franko.triangulation.SaveStructures.BaseFunctionValues;
import ua.lviv.franko.triangulation.SaveStructures.Node;
import ua.lviv.franko.triangulation.SaveStructures.NodePair;
import ua.lviv.franko.triangulation.SaveStructures.TriangleL;
import ua.lviv.franko.triangulation.SaveStructures.ValuesOnTriangle;
import ua.lviv.franko.visualisation.EFGResult2D;

/**
 * @author marko
 * 
 */
public class Solver {
	// general data
	protected int						quadraturePointsNumber;
	protected Gauss2D					gau;
	protected EquationParameters		parameters;
	protected HashMatrixLUStyle				mainMatrix;
	protected double[]					b;

	// data for background integration;
	protected ArrayList<TriangleL>		triangleList;
	protected ArrayList<Node>			nodeList;
	protected ArrayList<NodePair>		edges;
	// data for base functions
	protected ArrayList<Particle>		particles;
	protected ArrayList<BaseFunction>	baseFunctions;

	public Solver(int quadraturePointsNumber, ArrayList<Particle> particles, EquationParameters parameters, String folderDir, int PolinomialDegree, IWeightFunction2D weight) {
		// initialization
		this.quadraturePointsNumber = quadraturePointsNumber;
		this.gau = new Gauss2D(quadraturePointsNumber);
		this.particles = particles;
		this.parameters = parameters;

		// creation of base functions;
		this.createBaseFunctions(PolinomialDegree, weight);
		// continuation of initialization
		this.readData(folderDir);
		this.calculateValuesOnTriangles();
	}

	public EFGResult2D calculate() {
		for (int i = 0; i < b.length; i++) {
			for (int j = 0; j < b.length; j++) {
				System.out.print(mainMatrix.getValue(i, j) + "\t");
			}
			System.out.println();
		}
		for (int j = 0; j < b.length; j++) {
			System.out.println(b[j]);
		}
		for (int i = 0; i < b.length; i++) {
			BaseFunction bFunc = baseFunctions.get(i);
			// if
			// ((Math.abs(Math.abs(bFunc.getXi().getXi())-0.5)<0.0001)||(Math.abs(Math.abs(bFunc.getXi().getYi())-0.5)<0.0001))
			// {
			if ((Math.abs(Math.abs(bFunc.getXi().getXi()) - 0.5) < 0.0001)) {
				mainMatrix.setValue(i, i, 10000000000.0);
				b[i] = 0.0;
			}
		}
		LUFactorization lu = new LUFactorization(mainMatrix);
		return new EFGResult2D(lu.calculate(b), baseFunctions);
	}

	public void calculateValuesOnTriangles() {
		for (Iterator<TriangleL> iterator = triangleList.iterator(); iterator.hasNext();) {
			TriangleL triangle = iterator.next();
			// get nodes in which to calculate;
			ArrayList<Node> tmpNodes = new ArrayList<Node>();
			FuncWithJacobian fnc = new FuncWithJacobian(nodeList.get(triangle.A()), nodeList.get(triangle.B()), nodeList.get(triangle.C()), null);
			double by = 0;
			for (int i = 0; i < quadraturePointsNumber; i++) {
				by = 0.5 - gau.getX(i) / 2.0;
				for (int j = 0; j < quadraturePointsNumber; j++) {
					double[] xy = fnc.convert(0.5 + (0.5) * gau.getX(i), (by) / 2.0 + (by) / 2.0 * gau.getX(j));
					tmpNodes.add(new Node(xy[0], xy[1], false, j));
				}
			}
			// calculate function values and derivaties;
			ArrayList<BaseFunctionValues> bfvList = new ArrayList<BaseFunctionValues>();
			for (Iterator<BaseFunction> funcs = baseFunctions.iterator(); funcs.hasNext();) {
				BaseFunction function = funcs.next();
				BaseFunctionValues bfv = new BaseFunctionValues(tmpNodes.size());
				for (int i = 0; i < tmpNodes.size(); i++) {
					Node node = tmpNodes.get(i);
					double val = function.calculate(node.x, node.y);
					double valX = function.calculateDerivatyX(node.x, node.y);
					double valY = function.calculateDerivatyY(node.x, node.y);
					bfv.setValue(i, val);
					bfv.setDerivatyXValue(i, valX);
					bfv.setDerivatyYValue(i, valY);
				}
				bfvList.add(bfv);
			}
			triangle.setValues(new ValuesOnTriangle(tmpNodes, bfvList));
		}
	}

	public void generateMatrix() {
		mainMatrix = new HashMatrixLUStyle(RowType.ARRAY_TYPE, particles.size());
		b = new double[particles.size()];
		for (Iterator<TriangleL> iterator = triangleList.iterator(); iterator.hasNext();) {
			TriangleL triangle = iterator.next();
			for (int i = 0; i < b.length; i++) {
				// System.err.println(i);
				for (int j = 0; j < b.length; j++) {
					System.err.println(i + "\t" + j);
					double tmp = gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyXValue(), triangle.values.getValues().get(j).getFuncDerivatyXValue(),
							triangle.values.getNodes(), parameters.fKxx)
							+ gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyXValue(), triangle.values.getValues().get(j).getFuncDerivatyYValue(),
									triangle.values.getNodes(), parameters.fKxy)
							+ gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyYValue(), triangle.values.getValues().get(j).getFuncDerivatyYValue(),
									triangle.values.getNodes(), parameters.fKyy)
							+ gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyYValue(), triangle.values.getValues().get(j).getFuncDerivatyXValue(),
									triangle.values.getNodes(), parameters.fKyx);
					tmp += gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyXValue(), triangle.values.getValues().get(j).getFuncValue(),
							triangle.values.getNodes(), parameters.fPx)
							+ gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncDerivatyYValue(), triangle.values.getValues().get(j).getFuncValue(),
									triangle.values.getNodes(), parameters.fPy);
					tmp += gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncValue(), triangle.values.getValues().get(j).getFuncValue(),
							triangle.values.getNodes(), parameters.fQ);
					mainMatrix.setValue(i, j, mainMatrix.getValue(i, j) + tmp);
				}
				b[i] += gaussQuadratureIntegration(triangle.values.getValues().get(i).getFuncValue(), triangle.values.getNodes(), parameters.fF);
			}
		}
		for (int i = 0; i < b.length; i++) {
			for (int j = 0; j < b.length; j++) {
				System.out.print(mainMatrix.getValue(i, j) + "\t");
			}
			System.out.println();
		}
		for (int j = 0; j < b.length; j++) {
			System.out.println(b[j]);
		}

	}

	private void createBaseFunctions(int PolinomialDegree, IWeightFunction2D weight) {
		this.baseFunctions = new ArrayList<BaseFunction>();
		for (Iterator<Particle> iterator = particles.iterator(); iterator.hasNext();) {
			Particle particle = iterator.next();
			BaseFunction bf = new BaseFunction(PolinomialDegree, particle, weight);
			bf.setNodeList(particles);
			this.baseFunctions.add(bf);
		}
	}

	private double gaussQuadratureIntegration(double[] phi1, ArrayList<Node> nodes, IFunction2D coefFunc) {
		double res = 0.0;
		int pos = 0;
		double by = 0;
		for (int i = 0; i < quadraturePointsNumber; i++) {
			double temp = 0.0;
			by = 0.5 - gau.getX(i) / 2.0;
			for (int j = 0; j < quadraturePointsNumber; j++) {
				Node node = nodes.get(pos);
				temp += gau.getW(j) * phi1[pos] * coefFunc.calculate(node.x, node.y);
				pos++;
			}
			temp *= (by) / 2.0;
			res += temp * gau.getW(i);
		}
		res /= 2.0;
		return res;

	}

	private double gaussQuadratureIntegration(double[] phi1, double[] phi2, ArrayList<Node> nodes, IFunction2D coefFunc) {
		double res = 0.0;
		int pos = 0;
		double by = 0;
		for (int i = 0; i < quadraturePointsNumber; i++) {
			double temp = 0.0;
			by = 0.5 - gau.getX(i) / 2.0;
			for (int j = 0; j < quadraturePointsNumber; j++) {
				Node node = nodes.get(pos);
				temp += gau.getW(j) * phi1[pos] * phi2[pos] * coefFunc.calculate(node.x, node.y);
				pos++;
			}
			temp *= (by) / 2.0;
			res += temp * gau.getW(i);
		}
		res /= 2.0;
		return res;

	}

	private void readData(String folderDir) {
		triangleList = new ArrayList<TriangleL>();
		nodeList = new ArrayList<Node>();
		edges = new ArrayList<NodePair>();

		// read triangles;
		Scanner in = null;
		try {
			in = new Scanner(new FileInputStream(folderDir + "/file.1.ele"));
			int n = in.nextInt();
			in.nextInt();
			in.nextInt();
			for (int i = 0; i < n; i++) {
				in.next();
				triangleList.add(new TriangleL(Integer.parseInt(in.next()) - 1, Integer.parseInt(in.next()) - 1, Integer.parseInt(in.next()) - 1));
			}
		} catch (FileNotFoundException e) {
			e.printStackTrace();
		} finally {
			in.close();
		}
		// read nodes;
		try {
			in = new Scanner((new FileInputStream(folderDir + "/file.1.node")));
			int n = in.nextInt();
			in.nextInt();
			in.nextInt();
			in.nextInt();
			for (int i = 0; i < n; i++) {
				in.next();
				nodeList.add(new Node((Double.parseDouble(in.next())), (Double.parseDouble(in.next())), (isTrue(Integer.parseInt(in.next()))), i));
				// System.out.println(NodeList.get(i).x+"\t"+NodeList.get(i).y);
			}
		} catch (FileNotFoundException e) {
			e.printStackTrace();
		} finally {
			in.close();
		}
		// read edges;
		try {
			in = new Scanner((new FileInputStream(folderDir + "/file.1.poly")));
			in.nextInt();
			in.nextInt();
			in.nextInt();
			in.nextInt();
			int n = in.nextInt();
			in.nextInt();
			for (int i = 0; i < n; i++) {
				in.nextInt();
				int pos1 = in.nextInt() - 1, pos2 = in.nextInt() - 1;
				edges.add(new NodePair(nodeList.get(pos1), nodeList.get(pos2)));
				in.nextInt();
			}
		} catch (FileNotFoundException e) {
			e.printStackTrace();
		} finally {
			in.close();
		}

	}

	/**
	 * @param args
	 */
	public static void main(String[] args) {
		String FilePath = "C:/Users/marko/Documents/workspace/MyFEM/kvadrat/eps0";
		ArrayList<Particle> nodes = new ArrayList<Particle>();
		// параметри області:
		double ax = -0.5, bx = 0.5, ay = -0.5, by = 0.5;

		int n = 3, m = 3;

		double hx = (bx - ax) / (double) n, hy = (by - ay) / (double) m;

		// вузли:
		double rho = 2.1 * hx;
		for (int i = 0; i <= n; i++) {
			for (int j = 0; j <= m; j++) {
				Particle p = new Particle(ax + i * hx, ay + j * hy, rho);
				nodes.add(p);
			}
		}

		// equation parameters:
		EquationParameters params = new EquationParameters();
		params.setF(new ConstantFunction2D(10));
		params.setKxx(new ConstantFunction2D(1));
		params.setKyy(new ConstantFunction2D(1));
		params.setKxy(new ConstantFunction2D(0));
		params.setKyx(new ConstantFunction2D(0));
		params.setPx(new ConstantFunction2D(0));
		params.setPy(new ConstantFunction2D(0));
		params.setQ(new ConstantFunction2D(0));;

		Solver solver = new Solver(10, nodes, params, FilePath, 1, new SquareSupport(new CubicSpline()));
		solver.generateMatrix();
		solver.calculate().display(new Range(-0.5, 0.5), new Range(-0.5, 0.5));

	}

	private static boolean isTrue(int s) {
		if (s == 0) {
			return false;
		} else {
			return true;
		}
	}

}
